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ABSTRACT 

We present results of 2D3V particle-in-cell simulations of non-relativistic plasma collisions with 
absent or parallel large-scale magnetic field for parameters applicable to the conditions at young 
supernova remnants. We study the collision of plasma slabs of different density, leading to two different 
shocks and a contact discontinuity. Electron dynamics play an important role in the development of the 
system. While non-relativistic shocks in both unmagnetized and magnetized plasmas can be mediated 
by Weibel-type instabilities, the efficiency of shock-formation processes is higher when a large-scale 
magnetic field is present. The electron distributions downstream of the forward and reverse shocks 
are generally isotropic, whereas that is not always the case for the ions. We do not see any significant 
evidence of pre-acceleration, neither in the electron population nor in the ion distribution. 
Subject headings: acceleration of particles, instabilities, ISM:supernova remnants, methods:numerical, 
plasmas, shock waves 



1. INTRODUCTION 

Collisionless shocks in space environments are gener- 
ally considered to be sites of efficient particle accelera- 
tion. One of the main unsolved problems is how effi- 
ciently individual particles can be fed into a Fermi-type 
acceleration process. There, the particles have a mean 
free path, A m f p , large enough that they see the shock as 
a sharp discontinuity and assume a power-law spectrum 
(Blandford k Eichlcr 1987). By definition, collisionless 
shocks have a thickness of a few A m f p of the thermal ions. 
In particular for electrons, pre-acceleration appears to be 
required for injection into Fermi-type shock acceleration. 

The shock transition region harbors small-scale electric 
and magnetic turbulence that dissipate the flow energy 
and may provide pre-acceleration (Bret 2009). Studies 
of nonrelativistic shocks indicate that the structure and 
dissipation mechanisms at these shocks vary considerably 
depending on the upstream plasma parameters and the 
angle between the magnetic field and the shock normal, 
$Bn- Quasi-perpendicular shocks, 6>Bn > 45°, have been 
studied for many years (e.g., Papadopoulos 1988; Cargill 
k Papadopoulos 1988; Lembcgc et al. 2004), and the 
basic structure of such shocks with Alfven Mach num- 
ber Ma < 10 is well understood (e.g., Burgess k Sc- 
holer 2007). However, the physics of quasi-perpendicular 
shocks in the limit of high Alfven Mach number, as ap- 
plicable to supernova remnants (SNRs), is still under de- 
bate, in particular the electron injection problem (Umeda 
ct al. 2009; Lembege ct al. 2009; Amano k Hoshino 2010; 
Kato k Takabe 2010; Riquelme k Spitkovsky 2011; Mat- 
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sumoto et al. 2012). 

In contrast, quasi-parallel shocks, 6>B n < 45°, 
have been harder to understand. Observations of 
magnetic-field- aligned ion beams and large-amplitude 
low-frequency magnetic fluctuations upstream of space- 
plasma shocks have incited interest in electromagnetic 
ion-beam instabilities (Gary 1984) as sources of the 
necessary dissipation at the shock, particle injection, 
and subsequent scattering into a Fermi-type acceleration 
mechanism (e.g., Quest 1988). Studies of the shock tran- 
sition revealed that the shock ramp, where the principal 
density jump occurs, is not steady, but instead cyclicly 
re-forms such that the shock on average maintains its cor- 
rect speed (Winskc ct al. 1990; Krauss-Varban k Omidi 
1991). These effects have often been studied with one- 
dimensional hybrid simulations (Sugiyama 2011) that did 
not resolve electron and short-scale ion spatio-temporal 
scales. Short particle-in-cell (PIC) simulations of mildly 
rclativistic quasiparallel shocks were performed by Dieck- 
mann et al. (2010) and Murphy et al. (2010). In quasipar- 
allel conditions, Whistler and Alfven waves are excited 
that in their nonlinear evolution grow into Short Large- 
Amplitude Magnetic- field Structures (SLAMS), which 
are supposed to be efficient electron accelerators (Classen 
k Mann 1997; Stasiewicz k Ekeberg 2008; McClements 
k Fletcher 2009). 

Kato k Takabe (2008) recently reported the formation 
of nonrelativistic shocks in unmagnetized electron-ion 
plasmas via a Weibel-type short-wave instability (Wcibcl 
1959) for relative flow speeds between 0.2c and 0.75c (c 
is the speed of light). The relativistic two-stream Weibel 
(filamentation) instability (Medvedev k Loeb 1999) is 
known to mediate relativistic-shock formation in col- 
lisionless unmagnetized or weakly-magnetized electron- 
positron (e.g., Kazimura et al. 1998a; Silva et al. 2003; 
Nishikawa ct al. 2003; Chang et al. 2008) and electron- 
ion (e.g., Kazimura ct al. 1998b; Frcdcrikscn ct al. 2004; 
Spitkovsky 2008) plasmas. It leads to the development 
of current filaments and accompanying toroidal magnetic 
fields. The results by Kato k Takabe (2008) (see also 
Stroman et al. 2012) demonstrate that a pre-existing 
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magnetic field is not necessary for the existence of shock- 
mediating instabilities in non-relativistic flows. They 
also show that electron dynamics is relevant for shock 
formation in unmagnetized plasma, which calls for a 
fully kinetic (PIC) modeling in addition to hybrid simu- 
lations that can cover a longer time period and have been 
successfully employed for oblique nonrelativistic shocks 
(e.g., Gargate & Spitkovsky 2012). 

We have embarked on very-large-scale kinetic plasma 
simulations of the formation of unmagnetized or weakly- 
magnetized strictly parallel (#B n = 0°) nonrelativistic 
shocks, their long-term evolution, and particle accelera- 
tion. To address the question whether or not the shock 
structure and the particle spectra ever reach a steady 
state, we follow the evolution of the system longer than 
was done in Kato & Takabe (2008). We are interested 
in asymmetric flows, i.e., the collision of plasma slabs 
of different density, leading to two different shocks and 
a contact discontinuity (CD). The Mach number is set 
high enough that electrostatic shocks will not be formed. 
We want to see whether a parallel magnetic field ren- 
ders shock formation more efficient. In the simulations 
with magnetized plasma, the magnetic-field strength is 
assumed sufficiently low for the formation of the forward 
and reverse shocks with large Alfvenic Mach number, 
Ma ^> 1. This, together with the assumption of sub- 
relativistic plasma collision speed, renders our simula- 
tions applicable to the conditions at young SNRs. 

In contrast to a conducting-wall reflection method to 
set up the shock, used in Kato & Takabe (2008), our 
asymmetric slab-collision model is physically more accu- 
rate as it avoids assuming the existence of an infinitely 
sharp CD (in fact a grid-cell wide) and allows us to 
properly investigate the dynamics of the forward and re- 
verse shocks for different stream-counterstream parame- 
ters and magnetic field configurations. Beam-plasma 
systems are susceptible to various instabilities: some 
electrostatic (e.g., two-stream or Buneman modes), other 
quasi-electrostatic (e.g., Bret et al. 2010), and still others 
electromagnetic (e.g., filamentation) . The unstable spec- 
trum is thus at least two-dimensional. Which of these 
modes would grow fastest highly depends on the sys- 
tem parameters (Bret 2009). This clearly demonstrates 
a need for studies using setup like ours, since the most 
unstable modes excited in various settings may result in 
completely different shock-generating plasma waves. 

Besides demonstrating the efficiency of the shock- 
formation processes in the presence of the mean mag- 
netic field, a secondary advantage of restricting to a par- 
allel or absent large-scale magnetic field is that simula- 
tions can be set up cleaner. PIC simulations are sen- 
sitive to electromagnetic transients that arise from the 
initial setup and propagate through the simulation box. 
In shock simulations they can induce scattering turbu- 
lence and electron heating in the upstream region of the 
shock. Likewise, a perpendicular component of the ho- 
mogeneous magnetic field implies that in the collision 
region (or reflecting plane, if that is used to set up the 
shock) the curl of either the magnetic field or the mo- 
tional electric field is large, effectively acting as an an- 
tenna and emitting an electromagnetic pulse. 

The simulation model and setup are described in Sec- 
tion 2. The results of the simulations are presented in 
Section 3. We conclude with a summary and discussion 



in Section 4. 

2. SIMULATION SETUP 

In our simulations, two electron-ion plasma streams 
of different densities collide with each other at relative 
speed w rc i = 0.38 c. The density ratio between the dense 
and the tenuous plasma slab is 10. The simulation frame 
of reference is the center-of-momentum frame of the two 
plasmas, in which the dense (left) plasma stream moves 
to the right in -fx-direction with velocity vl = 0.0354 cx 
and the tenuous (right) stream moves to the left in — re- 
direction at vr = —0.354 ex. For simplicity, we will refer 
to the electrons and ions of the dense plasma stream as 
dense electrons and dense ions, and as tenuous electrons 
and tenuous ions for the dilute plasma. The plasmas 
are initially very cold; in the respective rest frames of 
the streams, the electrons have a thermal distribution 
with v e> th = y/kT/m e = 0.002 c and are in thermal equi- 
librium with the ions. We consider both unmagnetized 
and magnetized plasma conditions. In the latter case, 
the plasma flow is aligned with the homogeneous mag- 
netic field, i?o,:c, and thus we restrict our study to the 
formation of parallel shocks. The mean magnetic-field 
strength is given by the ratio of the electron cyclotron fre- 
quency Q e — eBo tX /m e to the electron plasma frequency 
of the dense-plasma component, w pe i = e 2 N e x / eom e 
(eo is the vacuum permittivity and N e< L is the density 
of dense electrons), which is f2 e /w P e,L = 0.04. Such a 
flow-aligned magnetic-field allows filamentation-like in- 
stabilities. As shown in, e.g., Bret et al. (2006), for 
the parameters of the initial plasma slab collision con- 
sidered here and taking into account counterstreaming 
electron beams only, a plasma magnetization larger than 
£l e /ujp et L ~ /3 ro i\/0T7rei ~ 0.125 would be required to 
suppress purely transverse filamentation modes, whereas 
a complete stabilization of oblique modes would never 
be efficient. With the assumed plasma magnetization of 
one-third of the suppression field, the filamentation insta- 
bility growth rate is reduced by 10% only. The simulation 
setup is thus designed to be applicable to young SNRs, in 
which dense supernova ejecta propagating with nonrela- 
tivistic velocity collide with a dilute weakly-magnetized 
interstellar medium. 

On account of the slower rate of shock-forming instabil- 
ities in nonrelativistic flows, our large-scale simulations 
of nonrelativistic plasma collisions have been performed 
in two spatial dimensions but keeping three components 
of the particle velocities (a 2D3V model). The code 
used in this study is a 2D3V modified version of the 
relativistic electromagnetic particle code TRISTAN with 
MPI-based parallelization (Buneman 1993; Niemiec et al. 
2008). We have performed extensive test simulations to 
formulate a numerical model that best conserves energy 
and minimizes numerical self-heating to provide a cold 
plasma beam that is sufficiently stable against numer- 
ical grid-Cherenkov and other instabilities. The model 
takes advantage of second-order particle shapes and uses 
a second-order FDTD field-solver with a weak Friedman 
filter (Greenwood et al. 2004) to suppress small-scale 
noise. In comparison with the numerical model used in 
previous work (e.g., Niemiec et al. 2008; Stroman et al. 
2012), that employed first-order particle shapes and a 
computationally-expensive iterative algorithm for digital 
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TABLE 1 

Parameters of the Main Simulation Runs. 



Run 


mi/m e 




Grid 








(^) 




Ul 


50 





6122.5x37.8 


4111 


U2 


20 





5531.7x39.8 


3803 


Ml 


50 


0.04 


6122.5x37.8 


4111 


M2 


20 


0.04 


5531.7x39.8 


3803 



Note. — Parameters of the main simulation runs described in 
this paper. Listed are: the ion-electron mass ratio rrii/m e , plasma 
magnetization, maximum computational box size in units of the 
ion skin depth of the dense plasma (X s i = 57.2 for runs Ul and 
Ml with mi/m e = 50 and X a i = 36.2 for runs U2 and M2 with 
mi/m e = 20), and the run duration in units of oj^ 1 . 



filtering of electric currents, the current model allows to 
maintain a satisfactorily low noise level with fewer parti- 
cles per cell, thus providing an increased overall perfor- 
mance of the code (compare Fonseca et al. 2008). 

Table 1 summarizes the parameters of the main sim- 
ulation runs discussed in this study. We use 5 particles 
per cell for each of the four plasma species in the ac- 
tive grid and apply the splitting method for the tenuous 
plasma, that assigns statistical weights to each particle 
to maintain the desired density ratio. Thus each ion par- 
ticle can be initialized at the same location as the cor- 
responding electron, to ensure a vanishing initial charge 
density. The electron skin depth of the dense plasma 
component X se .L — c/<jj pe> L — 8.1A, where A is the grid- 
cell size. Because the characteristic length scales of the 
plasma ions should be well contained in the simulation 
box, we can not use a realistic proton-to-electron mass 
ratio. Here we base our discussion on the results ob- 
tained for a reduced ion-electron mass ratio mi/m e = 50 
(runs Ul and Ml in Table 1). A slight dependence of the 
simulation results on the assumed mass ratio has been 
observed in Kato & Takabc (2008). To ascertain the 
effects of this parameter on the long-time evolution of 
the system we performed large-scale simulations also for 
mi/m e — 20 (runs U2 and M2); the results of these simu- 
lations are summarized in Appendix A. For runs Ul and 
Ml, referred here as to the unmagnetized and magnetized 
simulation, respectively, the ion skin depth of the dense 
plasma is X S i = X S i,L = c/u> P ix — 57. 2A. Hereafter, 
X s i will be used as the unit of length and we take the 
inverted ion plasma frequency of the dense component, 
uj^ 1 = u~^ L , as our unit of time. The time step we use 
is St = 0.062u; p ^ L = 0.0087a;-/. 

The simulations have been performed on a spatial grid 
in the x — y plane with radiating boundary conditions 
in x and periodic boundary conditions in y. The trans- 
verse size of the computational grid used in runs Ul and 
Ml is L y = 37.8 X s i and the system evolution has been 
followed for 4111a;- 1 , i.e., ~ 2.9 x lO 4 ^;-^, which in 

simulation Ml corresponds to 23.5 fi" 1 , where Oj is the 
ion cyclotron frequency. The simulations thus require a 
very large grid size in plasma flow direction, because the 
interpenetrating plasma streams occupy a continuously 
increasing volume. Our MPI-parallelized model parti- 
tions the simulation box into processor- domains of equal 
size, and the number of processor-domains in the trans- 



verse direction is constant. To limit the computational 
expense, in terms of both processor-usage and memory, 
an actual simulation is broken into several segments of 
different size in x, for which the number of processor- 
domains in flow direction is increased. Thus, to complete 
an simulation with, e.g., 200,000 time steps, succeed- 
ing runs with 1, 2, 4, and 6 216-processor domains each 
with 25,000 x 10 cells are performed on 216, 432, 864, 
and 1296 cores, respectively. Particles are continuously 
added to the system on both sides of the collision region 
with a so-called moving-injector. In this method, the in- 
jection layer moves away from the interaction region at 
a distance that preserves all the particles and fields gen- 
erated in the collision region, but at the same time does 
not allow a beam to travel a long distance without in- 
teraction with counterstreaming particles, which further 
helps avoiding numerical grid-Cherenkov effects. 

Note that the parameters used in the main numeri- 
cal experiments were chosen following an extended scan 
of the parameter space (e.g., number of particle per cell, 
the electron skin depth, the magnetic field strength) per- 
formed with smaller-size simulations. The results of the 
two selected production runs, Ul and Ml, are thus rep- 
resentative of the physical effects observed in the unmag- 
netized and weakly-magnetized system. The validity of 
the two-dimensional simulations in capturing the essen- 
tial physics of parallel shock formation through the elec- 
tromagnetic instabilities has been demonstrated in nu- 
merous studies of relativistic and nonrelativistic shocks 
(e.g., Spitkovsky 2008; Niemiec et al. 2008; Kato & Tak- 
abe 2008). 

The collision of the two plasma streams, as studied 
in this work, should in the nonlinear stage result in the 
formation of a system of forward and reverse shocks sep- 
arated by a CD. The shock structures observed in our 
simulations can thus be compared with the predictions 
of a hydrodynamic model. In analogy to the analysis pre- 
sented in Nishikawa et al. (2009), we can derive the CD 
speed, which in the simulation frame is /?cd = — 0.06x. 
The CD moves because the system is initially under mo- 
mentum balance but not ram-pressure balance. The ve- 
locities of the forward and reverse shock depend on the 
adiabatic index of the shocked plasma, which is different 
for unmagnetized and magnetized plasmas in the sim- 
ulations. In the unmagnetized case, particle motion is 
confined to the simulation plane and the nonrelativis- 
tic adiabatic index is T\j — 2. In the magnetized case, 
the magnetic field bends the particle trajectories out of 
the simulation plane. In effect, particles have three de- 
grees of freedom and the nonrelativistic adiabatic index 
is Tm = 5/3. In the following we use subscripts "U" 
and "M" to denote quantities measured in the unmagne- 
tized and magnetized run, respectively. Subscripts "R" 
and " L" refer to the rest frame of the tenuous and dense 
plasma, respectively. These subscripts are dropped for 
quantities measured in the simulation frame. 

For an unmagnetized plasma collision, the steady-state 
speed of the forward and the reverse shock should ap- 
proach /3fR ; tj = 0.43 and /3 r L,u = —0.14, in the rest frame 
of the tenuous and dense beam, respectively. In the sim- 
ulation frame, the shock velocities are /3f.u = 0.09 and 
/3 r .u = —0.11. The compression ratios at each shock are 
predicted to be Afu,d/^Vfu,u = 2.9 and iV r Tj,d/-/V r u,u = 
3.01. In the magnetized case, the simulation-frame ve- 
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locities of the forward and the reverse shocks should be 
/3f.M = 0.04 and t< m = —0.09, and the shock compres- 
sion ratios AffM,d/-NfM,u = 3.86 and A^ rM ,d/iV r M,u = 4.02, 
respectively. The Alfven velocityin the dense plasma is 
wax = [Bl x /^{N^ L m e + N^ L mi)] 1 / 2 = 5.67-10~ 3 c. In 
the rest frame of the tenuous plasma, the Alfven speed 
v a,r = 1-85 • 10~ 2 c. In the upstream rest frames of the 
tenuous and dense plasma the velocities of the forward 
and reverse shocks read /3fR,M = 0.39 and /3 r L,M = —0.13, 
respectively. Thus, the Alfvenic Mach number of the for- 
ward shock Mm = c(3{r,m/ v a,R ~ 21.1 and that of the 
reverse shock Mai — c/3 r L,M/«A,L — 22.5. 7 

3. SIMULATION RESULTS 

In what follows, we present the results of our two main 
numerical experiments for unmagnetized plasma condi- 
tions (run Ul; Sec. 3.1) and for magnetized plasma 
(run Ml; Sec. 3.2), both of which use an ion-electron 
mass ratio rrii/m e — 50. The results of two other large- 
scale runs, U2 and M2, that assume a lower mass ratio 
rrii/m e = 20, are described in Appendix A. The analysis 
of simulation data is also supported with additional test 
runs. 

To understand the nature of various unstable wave 
modes that appear in the systems, we also perform a 
linear analysis. Given the complexity of the situation 
and the need to investigate the full spectrum of unstable 
waves, our linear calculations are restricted to the cold 
plasma regime. We thus write the cold fluid equations 
for each of the plasma components and derive the di- 
electric tensor T(k, co) adapting the Mathematica Note- 
book described in Bret (2007). For both unmagnetized 
plasma and a flow-aligned magnetic field, we can with- 
out loss of generality restrict the wavevector to the x — y 
plane, k = (k x ,k y ). The dispersion relation then reads 
detT(k, uj) = 0, from which we calculate the growth 
rates of unstable wave modes. 

3.1. Unmagnetized Plasma 
3.1.1. Early-Stage Evolution 

In the first stage of the plasma collision, the interpene- 
tration of dense and tenuous streams quickly leads to cur- 
rent filamentation and the generation of transverse (out- 
of-planc in our two-dimensional set-up) magnetic field in 
the collision region. The observed turbulent modes result 
from Weibel-like instabilities, and their spectrum is dom- 
inated by modes with wave vectors oriented obliquely to 
the plasma flow (Bret ct al. 2005). This is in agreement 
with our linear analysis in Stroman ct al. (2012), that 
assumed similar stream-counterstream relative velocity 
and density ratio. The amplified magnetic field influ- 
ences the motion of particles. The ion beams slightly 
decelerate in bulk and develop a population of slow par- 
ticles. At the same time the dilute ion beam, in which 
filamentation is strongest, is heated. The bulk kinetic en- 
ergy released by the decelerating ions is converted to elec- 
tron heating in the electrostatic fields that accompany 
the ion filaments (e.g., Hedcdal ct al. 2004; Medvcdev 
2006; Nishikawa et al. 2006; Spitkovsky 2008). A fraction 
of the electrons of both streams is reflected away from 

7 A customary stated, the upstream flow Mach numbers read 
M A i = ?WfA,R - 19.1 and M Ar = ^l/^a.l - 6.3. 



the collision region. The initially reflected electrons of 
the tenuous plasma constitute a distinct hot beam that 
counterstreams against the incoming particles of the ten- 
uous stream. 



t=437u p; ~' 




3000 3050 3100 3150 

Fig. 1. — Structure of the plasma-collision region at time t = 
437 u)~- for the unmagnetized run Ul. Displayed are the profiles 
of (a) the average particle-number density normalized to the far- 
upstream density of the dense plasma (red lines: dense plasma, 
black lines: tenuous plasma; solid lines: ions, dotted lines: elec- 
trons), (b) profiles of the average magnetic (solid line) and electric 
(dotted line) energy density in simulation units, (c) the amplitude 
of the magnetic field B z (in sign-preserving logarithmic scale as 
sgn(ZJ 2 ) (3 + log [max(10 — 3 , |_B 2 1 )] )) , the density of dense electrons 
(d), dense ions (e), tenuous electrons (f), and tenuous ions (g), all 
normalized to the far-upstream density of the dense plasma, and 
the longitudinal phase-space distribution of electrons (h) and ions 
«■ 

By the time t ~ 400 w~ 1 , the magnetic turbulence has 
reached sufficient strength that fresh incoming upstream 
electrons can no longer cross the collision region. A con- 
tact discontinuity is formed that separates the electrons 
of both plasma slabs which slowly accumulate at this 
barrier. The structure of the collision region at this early 
evolution stage is shown in Figure 1. The CD is located 
around the initial collision boundary at x/\ S i ss 3120. As 
shown by the average-density profiles (Fig. la) and par- 
ticle density distributions (Figs, ld-lg), the CD creates 
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a barrier for less energetic dense ions. Only a few of these 
particles penetrate the tenuous-plasma region where they 
help maintaining charge neutrality. Since the number 
density of the dense ions in the tenuous-plasma region is 
comparable to that of the tenuous ions, the former play 
an important dynamical role in the system. The tenuous 
ions are only weakly perturbed upon crossing the CD. 
The fast tenuous-ion beam propagates into the extended 
region to the left of the CD. Its net charge is compen- 
sated in part by electrons of the dense beam that have 
been reflected during the initial plasma encounter and in 
part by hot tenuous electrons that leak across the CD. 
Very few dense electrons cross the CD into the tenuous- 
plasma region. In this zone, the charge balance with the 
dense ions streaming against the dilute beam is provided 
by reflected tenuous electrons. 

The structure of magnetic turbulence shows different 
characteristics on the two sides of the CD (Fig. lc). 
In the dense-plasma region, the tenuous fast beam con- 
stitutes a perturbation to the system. Ion-ion coun- 
terstreaming in this asymmetric beam/plasma system 
renders the unmagnetized oblique instability dominant 
in the spectrum of unstable modes. Oblique filamen- 
tary structures of wavelength A ~ 3.5A S i can be seen in 
the transverse component of the magnetic field and also 
through density modulations of the electrons and ions 
(Stroman et al. 2012). The density filaments are stronger 
in the tenuous plasma because they can be compensated 
by much weaker perturbations of the dense beam. Weak- 
beam conditions are also found at the leading edge of the 
dense ion beam propagating into the tenuous-plasma re- 
gion to the right of the CD. Here, the dilute beam is 
composed of the most energetic ions of the dense slab 
that have been accelerated in the collision region (Fig. 
li). The oblique mode is thus excited with wavelength 
A ~ 5A„i, comparable to the skindepth of the tenuous- 
ion component. However, the dense-ion beam shows a 
very strong density gradient in the tenuous-plasma re- 
gion. The conditions that determine the spectrum of 
unstable modes thus vary with distance from the CD 
from ones suitable for the oblique mode further upstream 
(at x/X S i w 3170) to those in the region closer to the 
CD (at x/\ si fa 3150) that allow quasi-perpendicular 
filamentation-like modes (kn pa 0) to grow fastest. In ef- 
fect, the magnetic structure formed at the dilute head 
of the dense-ion beam becomes strongly modified by a 
filamentation-like instability which grows faster than the 
oblique instability and therefore leads to more efficient 
magnetic- field amplification. As seen in Figure 1 , the in- 
stabilities spatially separate the ions of both slabs. Nar- 
row filaments in cold tenuous ions are surrounded by elec- 
trons. The volume between the filaments is depleted of 
tenuous ions and occupied by the ions of the dense beam. 
The filamentary structures in the electron components 
are more diffuse on account of the higher temperature 
of these particles. The turbulent magnetic structure ad- 
vected toward the CD is further amplified through merg- 
ers of current filaments. The observed characteristics of 
the magnetic-field generation in the tenuous-plasma re- 
gion are thus reminiscent of the Weibel-mediated rela- 
tivistic shocks in unmagnetized plasmas (e.g., Medvedev 
& Locb 1999; Silva et al. 2003; Nishikawa ct al. 2003; 
Frederiksen et al. 2004). 

Whereas the filamentation instability operating in the 



tenuous plasma can generate strong magnetic turbulence 
in this region, the oblique instability leads to only moder- 
ate magnetic-field amplification in the dense plasma left 
of the CD. As shown in Figure lb, the energy density 
of the magnetic field is more than an order of magni- 
tude larger in the region to the right of the CD, and this 
proportion is not significantly modified during the sub- 
sequent evolution of the system. The strong magnetic 
field generated in the CD region starts to efficiently de- 
flect and isotropize particles. The flow speed of the bulk 
plasma decreases, and larger numbers of ions are reflected 
back toward their respective upstream regions (Fig. li). 
The bulk kinetic energy of decelerating ions is converted 
into electron heating (Fig. lh). These processes of en- 
ergy dissipation should enforce the formation of a pair of 
shocks on the two sides of the CD. 

3.1.2. Late-Stage Evolution of the Forward-Shock 
Transition 

In the tenuous-plasma region, at t ~ 500 oT^- 1 the com- 
pression ratio reaches the level of ~ 3 predicted by the 
hydrodynamic jump conditions for a two-dimensional 
nonrelativistic plasma with adiabatic index T\j = 2. 
The compression of electrons is provided solely by the 
tenuous-electron component whereas the ion compression 
is initially dominated by the dense ions (see Fig. 2). The 
extended region between the unperturbed dilute plasma 
and the CD constitutes a forming forward shock. 

The structure of the forward-shock precursor (in the 
right half of the simulation box) is initially exclusively 
controlled by strongly nonlinear interactions between 
counterstreaming ion beams resulting in the formation 
and mergers of current filaments (see, e.g., Medvedev 
ct al. 2005; Polomarov ct al. 2008; Shvets ct al. 2009). 
The corresponding turbulence is not purely magnetic and 
transverse. The two ion beams pinch at a different rate, 
on account of their different density and temperature. 
Space charges thus arise that induce a weak in-plane 
electric-field, E yi which is not shown here (Tzoufras et 
al. 2006; Fiore et al. 2006). In the nonlinear phase, 
the magnetic-pressure gradient can also drive strong elec- 
trostatic fields (Dieckmann et al. 2009). The size of 
the magnetic precursor is thus determined by the range 
of the dense-ion beam in the tenuous plasma. Ahead 
of the filaments, a population of hot tenuous electrons, 
that have been reflected in the CD region, streams with 
mildly-relativistic velocity against the incoming beam. 
This excites a two-stream instability between reflected 
electrons and beam electrons that produces mainly lon- 
gitudinal (k x vr = 0) perturbations in E x and asso- 
ciated electron charge-density modulations with wave- 
length A ~ 5.5 \ S e,R, in agreement with our linear analy- 
sis. 8 The amplitude of electrostatic perturbations satu- 
rates at small level and their main effect is to heat-up the 
electrons of the incoming beam far upstream of the for- 

8 The linear analysis performed for conditions representative of 
the far-upstream forward-shock precursor also shows the existence 
of a weak filamentation mode that is in fact observed but does not 
provide a significant dynamical impact on the precursor. The third 
predicted mode — the Buneman mode resulting from the drift of 
electrons against the ions, that is needed for current balance with 
the reflected electrons, is an artefact of the zero-temperature limit 
assumed in the calculations and should dissipate quickly in the 
warm precursor plasma. 
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Fig. 2. — Structure of the forward-shock precursor at time t = 2405 w j . Shown are the profiles of (a) the average particle number 

density, (b) the average magnetic and electric energy density (see Fig. la-b), (c) amplitude of magnetic field B z (logarithmic scale; see 
Fig. lc), the density of tenuous electrons (d), tenuous ions (e), both normalized to the far-upstream density of the tenuous plasma, and 
the longitudinal phase-space distribution of electrons (f) and ions (g). 



ward shock. However, deeper in the forward-shock pre- 
cursor the electron-density perturbations become non- 
linearly enhanced through the approaching stream of 
the dense ions, and the fluctuations in E x are ampli- 
fied to amplitudes far exceeding those further upstream. 
This leads to a strongly nonlinear backreaction on the 
tenuous-plasma ions, whose density structure begins to 
resemble that of the electrons. Figure 2 presents the 
structure of the precursor region at time t = 2405 w" 1 
that is representative of the characteristics discussed 
above. The figure shows profiles of the average parti- 
cle densities (Fig. 2a) and electromagnetic energy den- 
sities (Fig. 2b), a two-dimensional distribution of the 
magnetic field B z (Fig. 2c), the densities of electrons 
and ions of the tenuous plasma (Figs. 2d and 2e), and 
the electron and ion phase-space density (Figs. 2f and 
2g). The electrons of the dense plasma are absent in 
this region. The density of the dense-ion beam is ho- 
mogeneous and shows a strong gradient (see red solid 
line in Fig. 2a and the beam with positive momentum 
in Fig. 2g). The far-upstream region is excluded from 
Figure 2, but the hot tenuous electrons counterstream- 
ing against the incoming beam can be clearly seen as a 
population with positive momenta beyond x/X si 3550 
in Figure 2f. The strong density fluctuations in the ions 
and electrons of the tenuous plasma can be seen in Fig- 
ures 2d and 2e for x/X S i > 3450. The longitudinal 
structures coexist with transverse density filaments in- 
dicating that the filamentation-like instability and the 



density modulation discussed above operate in parallel. 
Closer to the shock, the density modulations eventu- 
ally lead to the formation of plasma cavities. In this 
region, the density of the dense-ion beam is increased 
and also a population of reflected tenuous ions exists 
(for x/\ si < 3400; the density of the shock-reflected 
ions is more than a factor of 10 smaller than that of 
the dense ions). Such nonlinear structures are known to 
frequently occur in plasma physics and astrophysics ap- 
plications (e.g., Goldman 1984). As the parallel density 
structures are advected towards the shock and the den- 
sity of the dense-ion beam/shock-reflected ions increases, 
the modulations grow further in amplitude and eventu- 
ally begin to rapidly expand radially. The expansion and 
subsequent merging of the adjacent cavities causes a com- 
pression of the plasma between the voids and enhances 
electric currents in the compression regions that lead to 
additional magnetic-field amplification. 

The filamentary structure of the magnetic field thus 
becomes modified in the shock precursor. This can be 
observed in Figure 2c and in the mean magnetic-energy- 
density profile in Figure 2b (solid line) in which the ef- 
fects discussed are revealed as localized magnetic-field 
amplification. Note, that although the interaction be- 
tween the dense-ion beam and the charge-density waves 
turns highly nonlinear, it is the filamentation instability 
that remains the dominant source of magnetic-field gen- 
eration. As is evident in Figure 2b, the average magnetic- 
energy density in the region of nonlinear filamentation 
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{x/Xsi 3280) is about an order of magnitude larger 
that the energy density associated with plasma compres- 
sions further upstream. Nevertheless, the role of density 
modulations is not negligible. It is also important to note 
that the volume filled with weak electrostatic turbulence 
ahead of the magnetic-filamentation zone constitutes an 
inherent part of the shock precursor and as such must be 
included in realistic PIC modeling. The electrostatic pre- 
cursor is not a transient effect in nonrelativistic shocks, 
because there will always be a population of electrons re- 
flected upstream with relativistic velocities that outrun 
the bulk of the shock-reflected ions. 

In our symmetric unmagnetized test run (equal density 
in the colliding plasma slabs) we see similar effects in the 
shock precursor, but the cavities represent weaker den- 
sity fluctuations and do not grow nonlinearly as in the 
asymmetric simulation Ul. Consequently, magnetic-field 
filaments are only weakly modified by plasma compres- 
sion. The cavities that form upstream, dissipate closer 
to the shock and give place to density filaments. 



t=4i 1 iu p r' 




Fig. 3. — Structure of the collision region at the end of the sim- 
ulation Ul at time t = 4111 a;" 1 . Shown are the profiles of (a) the 
average particle-number density (a logarithmic scale is used), and 
(b) the average magnetic and electric energy density (see Fig. la- 
b), and the longitudinal phase-space distributions separately for all 
particle species: electrons of the dense (c) and tenuous (d) plasma 
and ions of the dense (e) and tenuous (f) component. In Fig. 3a, 
the horizontal dashed lines mark the hydrodynamic compression 
level of TVfu d/^fu.u = 2.9 for the forward shock (lower line) and 
^rU d/N r \j u = 3.01 for the reverse shock (upper line). 

Figure 3 shows the structure of the collision region at 
the end of the simulation at time t = 4111 wZ, . The 



profiles of the average particle densities (Fig. 3a) and 
electromagnetic-energy densities (Fig. 3b) are presented 
together with the phase-space distribution separately for 
all particle species (Figs. 3c-3f). The plasma density 
and electromagnetic-field structure (not shown) are anal- 
ogous to that at the earlier stage as described above 
(compare Figs. 2a-2b with Figs. 3a-3b). The CD is 
located at x/X S i « 3050. By the end of the run, its 
speed approaches the predicted /?cd = —0.06 in the sim- 
ulation frame. The forward shock region extends to right 
of the CD. The shock compression is in agreement with a 
hydrodynamic value of -/Vfu^/^Vfu.u = 2-9 (lower dashed 
line in Fig. 3a), but the shock front is poorly defined. 

The plasma filaments that grow in size when approach- 
ing the forward shock should disintegrate to form a 
nearly homogenuous downstream region. However, this 
stage is not yet observed by the end of the simulation. 
Although the downstream tenuous-ion density distribu- 
tion is strongly turbulent, it is evident from Figure 3f 
that the dilute-ion beam is not decoupled from the dense 
plasma and, despite being hot, continues to carry a sub- 
stantial bulk- flow energy across the CD. Consequently, 
the forward-shock compression ratio carries a substantial 
contribution from the dense ions. Note in Figure 3a, that 
the density of dense ions roughly equal to the tenuous- 
ion density is required to match the electron compression, 
that is provided solely by tenuous electrons. 
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Fig. 4. — Kinetic-energy spectra of particles downstream of the 
forward shock at x/X s i PS 3120 and time t = 4111 uj~^ for the 
unmagnetized run Ul, in the CD rest frame. In simulation units 
^kin = 0.25 (7 — 1) for electrons and -Ekin = 12.5 (7 — 1) for ions. 
Normalized spectra are shown for tenuous electrons (black solid 
line), tenuous ions (red solid line), and the dense ions (red dashed 
line) . 

As shown in Figure 4, that displays particle distribu- 
tions in the forward shock downstream, the tenuous-ion 
spectrum strongly deviates from a thermal distribution. 
The electrons of the tenuous plasma, accelerated at the 
expense of the ion bulk-kinetic energy, are thermalized. 
Since they are forced to follow the ion spatial distribu- 
tion to maintain charge neutrality, their density is neither 
uniform nor steady in the downstream region. We con- 
clude that the forward shock has not yet fully developed. 
Although in asymmetric unmagnetized-plasma collisions 
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kinetic instabilities generate ample magnetic turbulence 
in the tenuous-plasma region, the magnetic field is still 
too weak to efficiently decouple the dilute-ion beam from 
the dense slab and to force the ion-beam reflection needed 
to maintain the shock. 

The contribution of the reflected tenuous ions to the to- 
tal density of ions counterstreaming against the incoming 
tenuous beam and to the plasma compression at the CD 
continuously increases during the simulation. However, 
downstream of the forward shock the density of tenuous 
ions begins to exceed that of dense ions only at a very 
late stage. Reflection of the tenuous ions is most effec- 
tive downstream and in the foreshock region with strong 
current filaments (see region with 3100 < x/\ si < 3300 
in Fig. 3), but takes place also in the ion precursor fur- 
ther upstream where plasma cavities nonlinearly grow 
and expand (3300 < x/X S i < 3700). This can also be ob- 
served in Figure 2g. One can note that the locations at 
which the ions are decelerated and deflected well coincide 
with the density compressions caused by the expanding 
cavities (compare, e.g., x/X S i sa 3340 and 3390 in Figs. 
2d-e and 2g). As the electron phase-space plots reveal 
(Figs. 2f and 3d), the ion deceleration at these locations 
is accompanied by localized electron heating. 

As noted above, in the simulation frame the system 
of forward and reverse shocks separated by a CD moves 
toward the incoming dense-plasma beam with velocity 
/?cd = —0.06. This speed defines the CD frame or 
downstream rest frame, in which particle spectra are 
calculated. Particle spectra in kinetic energy, -Ekin = 
(7 — l)mi c 2 (mi — m e ,mi), are presented in Figure 4. 
They are based on all particles located downstream of 
the forward shock in a slice of width ~ 20 A s j, centered 
at x/X S i ~ 3120. At this location, the contribution of 
dense electrons (not shown) to the total number density 
and kinetic energy density is negligible. The distribution 
of tenuous electrons is isotropic and quasi-thermal with 
no evidence of a supra-thermal tail. Thus the efficiency 
of any electron-injection process that may have occurred 
in the forward-shock region in the unmagnetized plasma 
is very small (compare Sec. 3.2.5). The electrons ac- 
count for 37% of the total particle kinetic energy, dense 
ions contribute another 38%, and dilute ions carry the re- 
maining 25%. The ion distributions are very anisotropic, 
which impacts the energy spectra (compare Figs. 3e-3f). 
In the spectrum of tenuous ions, the decelerated beam 
largely contributes to a peak around -Ekin ~ 0.4; heated 
beam particles and particles returning upstream build 
the high-energy tail. The dense ions form a separate 
colder distribution that is not mixed in phase-space with 
the dilute component. The high-energy tail in the spec- 
trum is formed by particles that are reflected back to the 
dense plasma region. These features again demonstrate 
that the forward shock is yet in the formation phase and 
the processes of energy exchange/equilibration between 
the electron and ion components are still ongoing. 

3.1.3. Reverse-Shock Transition 

An important feature to be noted from Figure 3 is that 
the forward shock is a region in which the ions of the 
dense beam are strongly deflected and accelerated to- 
wards negative momenta, up to the initial momentum 
of the tenuous ions (see Fig. 4). In consequence, these 
particles are able to cross the CD back to the reverse- 



shock region to form there a population of hot and ener- 
getic ions counterstreaming against the incoming dense- 
plasma beam. This component of the dense ions can be 
clearly seen in the phase-space distribution displayed in 
Figure 3e. The counterstreaming dense-ion beam is gen- 
erated during the entire evolution of the system, slowly 
building the plasma compression associated with the re- 
verse shock left of the CD. 

By the time the simulation ends, the reverse-shock 
transition has not yet reached the compression ratio of a 
hydrodynamic shock (upper dashed line in Fig. 3a). Al- 
though the plasma density is homogeneous in the vicinity 
of the CD, particle distributions, including dense elec- 
trons, are strongly anisotropic. Upstream of the density 
compression, the oblique-mode instability amplifies the 
transverse magnetic field resulting in structures analo- 
gous to those observed earlier. The instability is gen- 
erated in the entire extended region that is populated 
with particles of the hot tenuous-ion beam. The mag- 
netic field is very weak and exerts no visible modula- 
tion to the plasma density. This is largely due to tem- 
perature effects. The tenuous-ion beam is accompanied 
by a cloud of very hot tenuous electrons that provide 
charge neutrality and also cause substantial heating of 
the dense electrons. The oblique-filamentary magnetic 
structures advected and compressed in the reverse shock 
are further amplified through interactions with energetic 
dense-ion particles reflected in the CD region, so that a 
weak magnetic precursor to the reverse shock forms. 

3.2. Magnetized Plasma 

The strength of the mean magnetic field assumed for 
the magnetized run Ml has been chosen to characterize 
the physical effects introduced to the system described 
in Section 3.1 by the presence of the homogeneous flow- 
aligned field component. The physical conditions repre- 
sent a weakly-magnetized plasma, in which strong shock 
waves should form. As checked with additional test runs, 
systems with lower magnetic-field strength show charac- 
teristics intermediate between those observed for run Ul 
and Ml. 

3.2.1. Early Evolution and Filamentation-like Instability 
in the Dense-Plasma Region 

Figure 5 shows the structure of the collision region in 
magnetized plasma at early time, t = 437 w~ 1 . As in the 
unmagnetized simulation (compare Fig. 1), the observed 
filamentation of plasma and the associated perpendicular 
magnetic turbulence result from Weibel-like instabilities. 
A direct comparison of Figures 1 and 5 reveals that sev- 
eral characteristics of the magnetized system are analo- 
gous to those described for the unmagnetized case. In 
particular, the magnetic-energy density in the tenuous- 
plasma region exceeds the turbulent-field energy reached 
in the dense plasma left of the CD. The similarities are 
also seen in the electromagnetic and density structure in 
the vicinity of the CD and in the tenuous-plasma region. 

Notably different from the unmagnetized case is the 
structure of the dense-plasma region, in which an almost 
purely perpendicular filamentation-like mode is produced 
at the head of the tenuous-ion beam counterstreaming 
against the incoming dense plasma. This mode appears 
in the out-of-plane magnetic field, B z , that is strongly 
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Fig. 5. — Structure of the collision region at time t = 437 w ]Ji 
the magnetized run Ml (see Fig. 1). Note that the mean magnetic 
field is subtracted in the magnetic energy-density profile in Fig. 
5b. The energy density associated with the mean field component 
is 7.9 X 10 — 3 in simulation units. 

amplified. The large-scale distributions of particle den- 
sities (Figs. 5d-5g) show that the density modulations 
associated with this mode first appear in the tenuous 
ions that form narrow filaments with large density con- 
trast (x/X si < 3020). Behind these structures, fila- 
mentary modulations develop in the densities of dense 
electrons and ions (3010 < x/\ si < 3040). These fil- 
aments merge and grow in size, but their density con- 
trast decreases further downstream of the tenuous-ion 
flow (3040 < x/\ S i < 3070), where magnetic-field fila- 
ments arc further amplified through mergers. 

The nature of the filamentation-like mode can be ver- 
ified with linear analysis. In Figure 6, we show the 
growth rate in reduced- wavevector space {Zn, Z±), Zi = 
kiVL/upe,L, for parameters that should represent the con- 
ditions at the head of the tenuous-ion beam. The cal- 
culations assume an initial density contrast of 0.1 be- 
tween the tenuous and the dense-ion beam and the cor- 
responding initial injection velocities. The charge and 
current of the dilute ion beam are to 90% balanced by 
electrons of the dense beam reflected during the initial 
plasma collision; the remaining portion is provided by 
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Fig. 6. — Linear growth rate (in units of the electron plasma fre- 
quency of the dense plasma, w~ e L ) of unstable modes as function 
of the flow-aligned and perpendicular wavenumber (fey and fcx), 
assuming parameters representing the conditions at the head of 
the tenuous-ion beam propagating in the dense plasma upstream 
of the reverse shock (run Ml). For a given k, the growth rate of 
the most unstable mode is plotted, and we present only the low-fc 
part of the spectrum that does not contain two-stream electrostatic 
modes. 

tenuous electrons leaking across the CD (see Fig. 5). 
The reflected dense electrons are assumed to move with 
the bulk of the tenuous beam to account for the cur- 
rent balance. Although the plasma beams are heated in 
the simulation, we perform our calculations in the zero- 
temperature approximation. However, in Figure 6 we 
present only the low-fc part of the spectrum, related to 
instabilities that involve the ions, and we do not show 
the two-stream electrostatic modes at high wavenumbcrs 
that grow strongly in a cold plasma but are suppressed 
at finite temperature. The perpendicular wavevector of 
the dominant unstable mode, Z±_ = 0.0138, is in per- 
fect agreement with Z± = 0.0125 — 0.0140 measured in 
the simulation (A ~ 2.25 — 2.5A S ;) and so is a nearly 
purely transverse character of the mode - the predicted 
Z±/Z\\ ~ 15- This suggests that the wavelength of the 
mode is not sensitive to temperature effects, although 
the growth rate appears to be reduced in warm plasma. 
We have further confirmed this with additional test sim- 
ulations that assumed higher initial temperatures of the 
colliding beams (see also Sec. 3.2.4). Our large-scale run 
M2, with rrii/m e = 20, together with another set of test 
runs, that incorporated larger ion-electron mass ratios 
of mi/m e — 100 and 200, reveal that a critical param- 
eter determining the strength of this unstable mode is 
ion inertia. The effect of the magnetic turbulence in the 
collision region on the motion of tenuous ions that pass 
through it is stronger for lower-mass particles and results 
in a warm beam for run M2 (see also Appendix A). In- 
stead, for larger mass ratios the physical conditions in 
the test simulations approximate a cold-plasma limit of 
the linear analysis. 

The presence of the filamentation-like mode signifi- 
cantly influences the dynamical evolution of the sys- 
tem. The nonlinear development of the instability in 
the dense-plasma region is presented in Figure 7 as a 
sequence of snapshots that show the B z magnetic-field 
component and density distributions of dense and ten- 
uous ions for tuj pi = 437, 655, 874, and 1093. The 
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Fig. 7. — Snapshots of the nonlinear evolution of the filamentation-like instability in the dense-plasma region upstream of the reverse 
shock at times (a) t = 437, (b) 655, (c) 874, and (d) 1093 a;~l; run Ml. From top to bottom, each panel shows: the amplitude of the 
magnetic field B z (in sign-preserving logarithmic scale; see Fig. lc), the density of dense ions, and the density of tenuous ions, both 
normalized to the far-upstream density of the dense plasma. 



distribution of dense electrons closely follows that of 
the dense ions. We also show in Figure 8 the aver- 
age profiles of particle densities and electromagnetic- 
energy densities, together with the phase-space distribu- 
tions across the whole collision region for tu> P i =1093 
and 1801. The filamentary modulations in the dense 
plasma, that were formed behind the transverse tenuous- 
ion density filaments during the linear phase of the in- 
stability (t = 437 w" 1 , Fig. 7a) turn in the later stage 
into dense-plasma cavities surrounded by regions of en- 
hanced plasma density. These structures form only in 
the dense plasma and can be explained as resulting from 
the magnetic bubbles that are created through localized 
accumulation of magnetic energy caused by merging of 
strong current filaments (see, e.g., Murphy ct al. 2010). 
The pressure gradient at the magnetic bubbles pushes 
the dense plasma out and causes the growth of the cav- 
ities. The subsequent radial expansion of the voids into 
each other results in further plasma compression and thus 
magnetic-field amplification. Similar in appearance to 
the dense plasma voids, narrow filaments in the tenuous 
ions, that initially emerge at the front of the beam, later 
merge, grow in size, and eventually disperse to develop 
a nearly homogeneous density distribution. The growth 
of the filaments can be traced through the magnetic-field 
structures in Figure 7. They demonstrate that the char- 
acter of the unstable mode evolves into an oblique mode 
that provides only a moderate magnetic-field amplifica- 
tion, compared to the level acquired through the initial 
filamentation-like mode (compare the average magnetic- 



field profiles in Fig. 5b for x/X S i < 3050 and in Fig. 8b 
for x/X S i < 2950). The tenuous ions remain unaffected 
by the pressure gradient at the magnetic bubbles due 
to their rapid motion. The density modulations in this 
plasma component seen around the cavities are instead 
the result of the nonlinear development of streaming in- 
stabilities that in this region become further enhanced 
by the approaching denser stream of tenuous electrons 
and dense ions reflected in the vicinity of the CD. This 
dense-ion component can be clearly seen in Figure 8d for 
3000 < x/X si < 3100. 

The formation of voids in the dense plasma greatly 
affects the subsequent evolution of the system. The cav- 
ities observed for 2970 < x/X si < 3050 at t = 1093 w^ 1 
(Fig. 7d) can be seen as modulations in the average den- 
sity of the dense plasma, that are associated with a broad 
peak in the magnetic-energy density (Fig. 8a-b) that is 
comparable to the magnetic-turbulence level generated in 
the tenuous-plasma region to the right of the CD. This 
strong field perturbs particle trajectories, thus enforcing 
the reflection of the interpenetrating beams into their re- 
spective upstream regions. One can observe this effect in 
Figures 8d-e. Localized particle reflection in the region 
of the cavities is independent of beam deceleration and 
reflection processes operating in the vicinity of the CD 
(3050 < x/X si < 3150). Due to the large strength of the 
magnetic field, which is continuously amplified by the 
compression between the expanding voids, particle re- 
flection is very efficient. One observes in the right panel 
of Figure 8 that by the time t = 1801 uj^ 1 returning ions 




Fig. 8. — Structure of the collision region in the magnetized simulation Ml, at time t = 1093<^> pi (left panel; a-e; compare Fig. 7d) 
and t = 1801 a;^ 1 (right panel; f-j). From top to bottom each panel shows: profiles of the average particle-number density, profiles of the 

average magnetic and electric energy density (see Fig. la-b), the longitudinal phase-space distribution for both electron species, and the 
phase space separately for the dense and tenuous plasma ions. In Figs. 8a and 8f, the horizontal dashed lines mark the hydrodynamic 
compression level of N{m d/^fM u = 3.86 for the forward shock (lower line) and iV r M d/^rM u = 4.02 for the reverse shock (upper line). 
The mean magnetic field with energy density of 7.9 X 10 — 3 is subtracted in the magnetic energy-density profiles in Figs. 8b and 8g. 



of both the dense and the tenuous plasma form almost 
homogeneous populations in the extended region around 
the CD. The region of cavities thus constitutes a barrier 
for the incoming beams that causes the accumulation of 
particles of each plasma component on their side of this 
barrier. The average-density profiles (Fig. 8f) and phase- 
space distributions (Fig. 8h-j) clearly illustrate this be- 
haviour. Note in particular the greatly increased contri- 
bution of the tenuous ions to the total density of ions that 
counterstream against the incoming tenuous beam right 
from the CD. Also, the dense-plasma region has become 
largely depleted of electrons of the tenuous component, 
that in addition have experienced strong heating in the 
region of cavities. 

The filamentation-like instability in the dense-plasma 
region thus provides a mechanism to efficiently decouple 
the colliding ion beams. Compared to the case of unmag- 
netized plasma, the presence of the mean magnetic field 
enables an efficient formation of the forward and reverse 
shocks, whose subsequent evolution is largely indepen- 
dent of the initial conditions. 

3.2.2. Late-Stage Structure of the Double-Shock System 

The structure of the collision region at the end of the 
simulation at time t = 41110;"/ is presented in Figure 9 



(to be compared with Fig. 3). One can clearly see a sys- 
tem of forward and reverse shocks separated by a CD. 
The shocks are visible as enhancement in the number 
density at x/X si ~ 3200 and x/\ si w 2830 for the for- 
ward and reverse shock, respectively. The compression 
ratio at each shock reaches ~ 4, in agreement with the 
hydrodynamical jump conditions for nonrelativistic gas 
with Tm = 5/3, that predict NfM,d/Nm.,u = 3.86 and 
-?VrM,d/A ir rM,u = 4-02 in the simulation frame for the for- 
ward and reverse shock, respectively (upper and lower 
dashed lines in Fig. 9a). The shock velocities measured 
in run Ml are also consistent with theoretical steady- 
state values of /^m = 0.04 and /3 r> M = —0.09 in the 
simulation frame. The CD moves with a speed close to 
the predicted (3cd — —0.06 and at the end of the run 
is located at x/X si ~ 2950. Note, that the location of 
the CD is shifted to the left compared to that in the 
unmagnetized run. Only in the magnetized case, the 
shock-formation processes allow for an efficient heating 
of the tenuous plasma, whose pressure is sufficiently high 
to drive bulk-plasma motion already at an early stage. 

As one can see in Figure 9a, the ion compression at the 
forward shock is almost entirely provided by the tenuous 
plasma (cf. unmagnetized simulation Ul). Although, as 
evidenced by the ion phase-space distributions in Figures 
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Fig. 9. — Structure of the collision region at the end of the simu- 
lation Ml at time t = 4111 w" 1 (see Fig. 3). The horizontal dashed 
lines in Fig. 9a mark the hydrodynamic compression level for the 
forward and reverse shock (see Fig. 8). The mean magnetic field 
with energy density of 7.9 X 10 — 3 is subtracted in the magnetic 
energy-density profile in Fig. 9b. 



9e and 9f, the contribution of the dense ions to the num- 
ber density of particles returning upstream is comparable 
to that of the tenuous ions, the latter represent a more 
directed and faster flow deep in the forward-shock precur- 
sor and also completely dominate the precursor dynamics 
further upstream (beyond x/X S i s» 3700). Thus, late in 
the simulation the forward shock begins to resemble a 
self-propagating structure that is exclusively maintained 
by shock-reflected particles, independent of the processes 
operating at and behind the CD. This is not yet the case 
in the dense-plasma region. Although the dynamics of 
the deep reverse-shock precursor (2750 < x/X si < 2820) 
is solely governed by the shock-reflected dense ions, the 
physics farther upstream is still determined by the inter- 
action of counterstreaming beams of dense and tenuous 
ions (see below). Note, however, that the tenuous-ion 
beam in the reverse-shock precursor, that has diluted, 
slowed-down, and heated, is at this stage decoupled from 
the main tenuous-plasma beam coming in from the right 
of the CD (Fig. 9f). 

The turbulent electromagnetic field in the collision re- 
gion is strongly amplified throughout the simulation to 
eventually reach a level about an order of magnitude 
larger than that achieved in unmagnetized conditions. 
One can see in Figure 9b that the turbulence is almost 
purely magnetic, and both the forward and the reverse 



shock develop strong magnetic precursors. The energy 
density of the magnetic field in the forward-shock region 
peaks at the shock front and then decays downstream 
(Kato & Takabc 2008). This turbulence is mostly gen- 
erated when the forward shock is self-propagating. The 
second maximum to the immediate right of the CD arises 
from the amplified magnetic turbulence produced via fil- 
amentation instabilities during the early evolution of the 
system, before the ion beams have become decoupled 
from each other. The broad peak in the magnetic field 
energy density to the left of the CD has largely evolved 
from the magnetic turbulence generated earlier in the 
region of cavities (Fig. 8). It also constitutes the shock- 
compressed magnetic turbulence produced upstream of 
the reverse shock after the cavities in the dense plasma 
have dissipated to form the downstream region of the 
reverse shock. 

3.2.3. Structure of the Forward-Shock Transition 

Figure 10 presents the structure of the forward-shock 

transition at time t = 3499 • This time was chosen to 

pi 

display the magnetic-field amplitude B z (Fig. 10a) and 
the density of tenuous ions (Fig. 10b) with correct aspect 
ratio; the characteristics discussed are the same as those 
observed at t = 4111 (Fig. 9). The tenuous-electron 
distribution (not shown) closely follows that of tenuous 
ions; they also provide charge balance to the low-density 
dense ions that occupy the volume between tenuous-ion 
filaments. As seen in Figure 10, the main source of tur- 
bulent magnetic field is the electromagnetic Weibel-type 
instability. Strong nonlinear effects of the interaction of 
returning particle beams with electrostatically-induced 
parallel density structures, reported in the forward-shock 
precursor for the unmagnetized case, are not observed. 
The two-stream instability between incoming and return- 
ing electron beams generates electrostatic perturbations 
far upstream that lead to nonlinear charge-density waves 
in the precursor (for x/X si > 3450 in Fig. 10b), but 
the interaction of the density modulations with the re- 
flected ions does not result in the formation of plasma 
cavities. The growth and expansion of parallel density 
fluctuations is suppressed by the mean magnetic field. 

The tenuous-plasma filaments, that pinch and merge 
in the forward-shock foot, break up at the shock inter- 
face and dissipate to form a nearly homogeneous density 
distribution downstream of the shock (x/\ s i < 3150 in 
Fig. 10b). The disintegration of the tenuous-ion beam 
downstream of the forward shock was not observed in the 
unmagnetized simulation Ul. Nevertheless, although the 
downstream tenuous-ion density is close to uniform, they 
are not fully thermalized (see below), which indicates 
that the forward-shock evolution has not fully reached a 
steady state by the end of the simulation. 

3.2.4. Structure of the Reverse-Shock Transition 

As discussed in Section 3.2.1, strong fluctuations in the 
density of the dense plasma and the magnetic field result- 
ing from a filamentation-like instability in the region of 
cavities began to form an obstacle for the incoming ten- 
uous plasma and eventually caused a depletion of these 
particles in the dense-beam region upstream of the dense- 
plasma voids. In addition, the tenuous ions at the head 
of the dilute-ion beam were heated and the filaments dis- 
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Fig. 10. — Structure of the forward-shock precursor in run Ml at time t = 3499 u) - . Shown is (a) the amplitude of the magnetic field 



B z , (in sign-preserving logarithmic scale; see Fig. 
dilute plasma. 



lc), and (b) the density of tenuous ions normalized to the far-upstream density of the 



persed. In effect, the character of the unstable magnetic 
mode has changed from almost purely perpendicular to 
slowly-growing oblique turbulence, with characteristics 
similar to those described for the reverse-shock precur- 
sor in the unmagnetized run Ul (Sec. 3.1.3). As in the 
unmagnetized case, the oblique magnetic mode is gener- 
ated wherever the tenuous-ion beam propagates in the 
upstream region of the reverse shock. However, the pres- 
ence of the mean magnetic field triggers an additional 
plasma mode in the system that slowly grows and even- 
tually dominates the magnetic-turbulence spectrum in 
the reverse-shock precursor. 
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Fig. 11. — Structure of the reverse-shock transition in run Ml at 
time t = 3061 wZj . Shown arc (a) profiles of the average particle- 
number density (see Fig. la) and the amplitude of magnetic-field 
components, B y (b), and B z (c). The horizontal dashed line in Fig. 
11a marks the hydrodynamic compression level of iV r M d/^iM,u = 
4.02 for the reverse shock. Note that the field components are 
in sign-preserving logarithmic units, e.g., for B z as sgn(B z ) (3 + 
log [ma X (10- 3 ,|B z |)]). 

Figure 1 1 presents the reverse-shock transition at time 
t = 3061 w" 1 , at which the initial oblique mode has al- 
ready attenuated. As seen in Figures lib and 11c, that 
show the By and B z magnetic-field components, respec- 
tively, the unstable mode represents an almost purely 
parallel magnetic wave. The wave is left-hand circu- 
larly polarized 9 and in the rest frame of the dense-plasma 
beam moves antiparallel to the mean magnetic field with 



phase velocity v pn ~ 7.9 x 10 



1.4^a,l- Therefore, 



in the frame of the reverse shock the wave is carried to- 
wards the shock with the plasma flow. The wavelength of 

9 For Bo = Bo,xX and k = k x x, a plane electromagnetic wave 
of the form Bi expi(k ■ x — bit), Bi = -Bi, y y + -Bi.z z is left-hand 
circularly polarized if -Bi, z = —iB\ tV . 



the mode is A ~ 15 X S i and its frequency is uj ~ 0.58 fV 
These characteristics do not significantly vary over the 
precursor region. 

To understand the nature of the observed parallel 
mode, a linear analysis has been performed. It is anal- 
ogous to the analysis for the filamentation-like mode in 
Section 3.2.1, but uses a modified set of parameters, for 
the physical conditions in the reverse-shock precursor 
have changed. The calculations thus assume a density 
contrast of 0.05 between the dilute and the dense-ion 
beam (compare the average density profiles in Fig. 11a) 
and a tenuous-ion beam velocity of v = —0.3 c, slightly 
slower than the initial injection velocity Vr. The charge 
and current of the tenuous ions is compensated mainly 
by returning dense electrons; electrons of the tenuous 
plasma represent only 5% of all returning electrons. We 
assume the cold-plasma approximation. The results are 
displayed in Figure 12, in which the two-stream electro- 
static modes at high wavenumbers are again not shown. 
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Fig. 12. — Linear growth rate of unstable wave modes for con- 
ditions representing the late-stage structure of the reverse-shock 
precursor in run Ml. (see Fig. 6). 

The area of strongest growth at Zii ~ 0.0019 and 
Z± s» 0.012 pertains to an oblique mode, analogous to 
the nearly transverse mode discussed in Section 3.2.1. 
The second maximum at Z± s» and Ztt s» 0.002, or 
Aji ~ 16A S i, can be identified with our parallel mode. 
Note, that the growth rates of both the oblique and 
the parallel modes are comparable to within 6% in the 
present setup, which was not the case for the analysis 
of the early-stage turbulence modes in the dense-plasma 
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region, in which a parallel mode can also be observed at 
Z\\ ss 0.0012 but with a growth rate 70% of that of the 
dominant unstable filamentation-like mode (Fig. 6). 

The nature of the parallel mode can be explained by 
analysis of the dielectric tensor T(k, w) for our system. 
For flow-aligned perturbations with k y = 0, the tensor T 
takes the form, 

/ b -c\ 
T= a , (1) 
V c b J 

where the tensor elements a, 6, c are reported in Ap- 
pendix B. As one can see, the dispersion equation dis- 
plays three branches, as detT(k, u>) = readily gives 

a(b + ic)(b - ic) = 0. (2) 

The equation a = defines two-stream like electronic 
modes at large wavevectors. As already noted, at 
the time-scale of the present simulations, these modes 
quickly saturate for warm electron beams. On the 
ion time-scale, the large-wavelength modes governing 
the system arise from the branches (b + ic) — and 
(b— ic) = 0. When analyzing a simpler setup formed by a 
dilute-ion beam passing through an electron-ion plasma, 
the so-called nonresonant (Bell) mode, as well the res- 
onant one, both come from (b — ic) =0 (Bell 2004; 
Reville 2006; Winske & Leroy 1984). In the present 
case, the growth-rate profile observed along the parallel 
axis in Figure 12 (and also in Fig. 6) can equally be 
traced back to the same (b — ic) = branch. We can 
therefore identify our parallel mode with the "Bell- like" 
mode. However, the linear analysis predicts a frequency 
to = 0.0012 uj pe x = 0.0085 LOpi^L — 1-5 f2< 9 that is a factor 
of 2.6 higher than observed in our simulation. This indi- 
cates that either the cold plasma approach is invalid or 
the simplified conditions assumed for the linear analysis 
do not properly represent the conditions in the simula- 
tions. 

To evaluate how well the simplified conditions corre- 
spond to the state of the plasma system in the simulated 
reverse-shock precursor and to estimate the effects of a 
finite temperature, two additional test runs have been 
performed. Both assume the conditions used in the lin- 
ear analysis and periodic boundary conditions. In the 
first run all plasma components are cold. We refer to 
this simulation as to the cold test-run. In the second test 
simulation, the plasma beams initially have temperatures 
as measured in the precursor. This run is referred to as 
the warm test-run. We found that both test runs well 
reproduce the characteristics of the initial oblique mode, 
which is the fastest-growing mode in the linear analysis. 
The nearly homogeneous density distributions, and the 
absence of strong density filaments, that are observed in 
the main simulation, are both well replicated in the warm 
test-run. In the latter and in the magnetized run Ml, the 
oblique mode is observed mainly through filamentation 
in the B z field component. The parallel mode appears 
at a later stage in both test runs. Whereas in the warm 
run it is the only magnetic mode, in the cold test-run 
the oblique-mode nonlinear filaments in B z are clearly 
imprinted onto the parallel magnetic- field fluctuations. 
These characteristics are due to thermal effects that lead 
to damping of the wave modes with perpendicular com- 



ponents, because the plasma pressure opposes pinching 
of the filaments. Consequently, the oblique mode even- 
tually dissipates in the warm plasma, and the parallel 
mode, whose linear growth rate is only slightly lower 
than that of the oblique mode, dominates the nonlinear 
evolution. 

The parallel mode observed in the test runs has a wave- 
length Aij ~ 16A S! ; {Z\\ ~ 0.002), is left-hand circularly 
polarized, and propagates antiparallcl to the mean mag- 
netic field. The measured growth rate is about half the 
analytically-derived growth rate, which is consistent with 
the influence of thermal effects. These results very well 
reproduce the properties of the parallel mode observed in 
the magnetized simulation. The measured frequency of 
the mode isww 1.1 fli in the cold test run and ui ~ 0.7 fli 
in the warm test-run. The results of the warm test run 
are therefore consistent with our main simulation, which 
indicates that the simplified conditions well represent the 
real conditions in the simulations and thus our identifi- 
cation of the parallel mode with the Bell-like mode is 
correct. The differences between the measured mode fre- 
quency and that predicted in the linear analysis likely 
result from thermal effects, that in the nonlinear stage 
modify the physical conditions, possibly introducing fac- 
tors not accounted for in the linear calculations. 

As one can note in Figure 11, the amplitude of the par- 
allel magnetic-field fluctuations considerably increases in 
the reverse-shock foot. The latter is the precursor re- 
gion in which a high-density population of shock-reflected 
dense ions is present at x > 2880 X S i. The amplifica- 
tion and compression of the parallel fluctuations thus 
provides a formation mechanism for the reverse-shock 
transition that is consistent with the hydrodynamic pic- 
ture. Note, that although the parallel magnetic mode 
excited in the precursor region results from the initial 
conditions, similar unstable modes should appear later 
on account of shock- reflected dense plasma. In fact, our 
linear calculations show that the appearance of the Bell- 
like mode critically depends on the presence of returning 
dense electrons in the shock-precursor plasma, and not 
on the tenuous-ion beam. Therefore, one may expect 
that the further evolution of the reverse-shock system 
will be similar. 

3.2.5. Particle Distributions and Electron Injection at the 
Forward Shock 

As in the case of unmagnetized plasma (run Ul), 
particle distributions are calculated in the rest frame 
of the downstream plasma, defined by the CD speed 
/?ce> — —0.06 in the simulation frame. Although in run 
Ml the downstream plasma is turbulent on large scales 
in the CD rest frame, these turbulent motions are slow 
and do not significantly distort the spectra. 

For run Ml we define the forward-shock downstream 
region as the area of roughly constant plasma compres- 
sion with ratio ~ 4, that is dominated by particles of the 
tenuous plasma whose bulk velocity is similar to that of 
the CD. At t = 4111 c^ 1 , the downstream region spans 
the range 2980 < x/X si < 3080 (see Fig. 9). Left of 
this location is the CD region in which both the dense 
and the tenuous plasma reside and where additional pro- 
cesses may operate that alter the particle spectra. As 
for run Ul, to calculate particle distributions we choose 
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a slice of width ~ 20X si , now centered at x/X S i ~ 3000. 
The slice is stationary with respect to the CD, and thus 
the shock front moves steadily away from the location at 
which spectra are measured. 
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Fig. 13. — Kinetic-energy spectra of tenuous electrons (black line) 
and tenuous ions (red line) in the forward-shock downstream region 
at x/\ s i R2 3000 and time t = 4111 u>~^ for the magnetized run Ml, 
in the CD rest frame. The spectra are normalized and expressed 
in simulation units, in which m e c 2 = 0.25. 

Figure 13 presents spectra of tenuous electrons and 
ions at time t — 4111 . Only the distributions of elec- 
trons appear approximately thermalized. Visible also in 
Figure 9, the distribution function of tenuous ions ex- 
hibits a plateau when integrated over any two momen- 
tum coordinates, e.g., f(p x ) = J dp y dp z f(p), and con- 
sequently the kinetic-energy spectrum strongly deviates 
from a thermal distribution. Isotropy is approximately 
maintained in the plateau section of the distribution, but 
not in the wings where the distribution function falls off. 
The form of the distribution suggests that the ions have 
reached stability against parallel electrostatic and elec- 
tromagnetic instabilities. As in the case of the unmag- 
netized system, in run Ml the total energy content in 
electrons is roughly the same as that in ions, in fact 2/3 
of it. 

Figure 14 details the time evolution of the downstream 
electron spectrum. At all times for which spectra are 
presented {tw pl = 2406, 2974, 3542, and 4111, shown 
with black, blue, green, and red line, respectively), the 
forward-shock compression is well established in the ex- 
tended region to the right from the CD, and the elec- 
trons are isotropic in the downstream frame (cf. run Ul). 
These spectra display a slow increase in average energy, 
possibly because the region, in which the distribution is 
calculated, is located closer to the shock front at earlier 
times. All spectra are slightly more peaked than a rel- 
ativistic Maxwellian, which we plot for comparison. It 
is interesting to note that renormalizing the spectra to 
dimensionless energy x — £ , /£'km,av, £kin,av being the 
average particle energy, indicates that the shape of the 
spectra evolves very weakly, and it is indeed mainly a 
continuous increase in -Ekin.av that is responsible for the 
spectral evolution. There is no significant indication of 
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Fig. 14. — Tenuous-electron kinetic-energy spectra in the 
forward-shock downstream region at ~ 50A si from the CD at times 
t = 2406 (black line), 2974 (blue line), 3542 (green line), and 
4111 (red line) for run Ml. All spectra are normalized and 

expressed in simulation units, in which m e c 2 = 0.25. The dotted 
line shows a relativistic Maxwellian for comparison. 

a supra-thermal tail in the electron spectrum. 
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Fig. 15. — Kinetic-energy spectra of dense electrons (black line) 
and dense ions (red line) in the reverse-shock downstream region 
at x/X s i ft: 2850 — 2890 and time t = 4111 uj~^ for the magnetized 
run Ml, in the CD rest frame. The spectra are normalized and 
expressed in simulation units, in which m e <? = 0.25. The dotted 
line shows a Maxwellian fit to the electron spectrum. 

Dense-plasma particles downstream of the reverse 
shock are sampled in a slice of width ~ 40 A si , centered 
at x/X S i ~ 2870 and their energy distributions in the CD 
rest frame at the end of run Ml are presented in Figure 
15. The hydrodynamic plasma compression is fully es- 
tablished in this narrow region at the very end of the 
simulation. The electrons and ions carry to within 10% 
the same kinetic energy. Only in the high-energy tail do 
we observe a modest anisotropy in the electron distribu- 
tion. A Maxwellian fits the electron spectrum reasonably 
well; a slight discrepancy at high energies coincides with 
anisotropy. The ions are anisotropic as well; a substantial 
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population of returning particles exists that builds a su- 
perficial high-energy tail (see Fig. 9), which is therefore 
not indicative of ion pre-acceleration. 

4. SUMMARY AND DISCUSSION 

We have performed 2D3V PIC simulations of non- 
relativistic plasma collisions with absent or parallel large- 
scale magnetic field. Contrary to earlier expectations, 
simulations by Kato & Takabc (2008) suggested that 
Weibel-type instabilities can trigger shock formation in 
an unmagnetized environment even for collisions with 
non-relativistic velocity, but fast enough to prohibit an 
electrostatic shock (Forslund & Shonk 1970; Forslund 
& Frcidbcrg 1971). A pre-existing magnetic field may 
therefore be not necessary for the formation of non- 
relativistic shocks. The question we addressed in this 
study is whether or not the shock structure and the par- 
ticle spectra reach a steady state, and to what degree a 
parallel magnetic field renders shock formation more effi- 
cient. A secondary goal was to investigate the efficacy of 
particle pre-acceleration, which appears needed for their 
injection into a diffusive shock-acceleration process. 

To address these questions we follow the evolution of 
the system longer than was done by Kato & Takabe 
(2008). Also, we are interested in asymmetric flows, i.e. 
the collision of plasma slabs of different density, leading 
to two different shocks and a contact discontinuity which 
is self-consistently modeled. Emphasis has been placed 
on setting up the simulations as clean as possible. 

The simulation parameters are chosen to be close to 
those of young supernova remnants, except that we use a 
reduced ion-to-electron mass ratio, m,/m e = 50. About 
50 years after the explosion of a type-la supernova, the 
density ratio between ejecta and ambient medium is 
about 10 as assumed here, if the density of the ambient 
gas is about 1 atom/cm 3 . The shock velocities are about 
0.05 c, a factor of a few lower than in our simulations. In 
the magnetized case, the Alfvenic Mach number of the 
forward shock is M^f = c/3^i m/ v a,k — 21.1 and that 
of the reverse shock Ma? = c/3 r L,M/vA,L — 22.5. Here, 
as forward shock we denote the shock propagating into 
the low-density (ambient) plasma; vice versa the reverse 
shock travels through the dense plasma, that represents 
the ejecta. 

Our results can be summarized as follows: 

• Both in the unmagnetized and the magnetized case, 
a double-shock structure builds on time scales of 
about a thousand ion plasma times, w" 1 . The 
colliding plasmas are initially very cold, so that 
the sonic Mach number in our simulations is high 
enough to prevent an electrostatic shock. We show 
that while non-relativistic shocks in both unmag- 
netized and magnetized plasmas can be mediated 
by Weibel-type instabilities, the efficiency of shock- 
formation processes is higher when a large-scale 
magnetic field is present. The amplitude of mag- 
netic turbulence generated in the collision region 
is in the magnetized plasma an order of magni- 
tude larger than in unmagnetized conditions. The 
typical magnetization of the ambient medium, 
into which young SNR shocks propagate, given by 
ct = MJ 2 , is in the range 10~ 6 < a < 10~ 3 . Our 
magnetized run with a « 2 x 1Q -3 thus probes an 



upper limit of the sigma parameter, and any real- 
istic formation scenario of parallel shocks in SNRs 
should show characteristics intermediate between 
the two cases studied here. 

• In the unmagnetized case, small-scale electrostatic 
and filamentation-like instabilities operate in par- 
allel and heat the electrons. Eventually, strong 
fluctuations arise in the density of electrons and 
ions that lead to the formation of the double-shock 
structure. The electron dynamics play an impor- 
tant role in the development of the system. Ion-ion 
or ion-electron streaming generally drives the tur- 
bulence, which is mainly magnetic. The exact type 
of an instability, however, depends on the location 
in the shock precursor and is generally different for 
forward and reverse shocks. Nonlinearities and the 
change of conditions during plasma advection to- 
ward the shocks render an identification with linear 
modes difficult. In the magnetized case, filamenta- 
tion in the dense plasma region is initially more 
pronounced, which we verify with linear analysis, 
but it quickly evolves into an oblique mode and 
nonlinear density fluctuations. Some of the stream- 
ing instabilities that lead to the formation of the re- 
verse shock are similar to non-resonant modes that 
have been discussed in the context of cosmic-ray in- 
duced magnetic-field amplification (e.g. Bell 2004; 
Niemiec et al. 2008; Riquelme & Spitkovsky 2009; 
Stroman et al. 2009; Niemiec et al. 2010; Stroman 
ct al. 2012). In both the magnetized and unmagne- 
tized systems, magnetic-field generation processes 
lead to stronger magnetic fields downstream of the 
forward shocks, rather than reverse shocks. We 
would thus expect synchrotron emission to origi- 
nate from the magnetized structure of the forward 
shock transition, which is in agreement with SNR 
observations. 

• The electron distributions downstream of the for- 
ward and reverse shocks are generally isotropic, 
whereas that is not always the case for the ions. 
The observed spectra suggest that the electrons 
have reached a statistical equilibrium and gener- 
ally carry nearly as much kinetic energy as do the 
ions, although the latter have not yet reached an 
equilibrium. 

• We do not see any significant evidence of pre- 
acceleration, neither in the electron population nor 
in the ion distribution. Note that the shock-surfing 
acceleration, that reportedly provides fast electron 
acceleration (Amano & Hoshino 2009), should not 
operate in situations with absent or strictly par- 
allel magnetic field, because v x B = 0. In our 
magnetized simulation, the Alfvenic Mach number 
is slightly below the threshold (at Ma > 25) for the 
electron-acceleration process described by Amano 
& Hoshino (2010). Dieckmann et al. (2010) and 
Murphy et al. (2010) report efficient electron ac- 
celeration to Lorentz factors of order 100 in PIC 
simulations of slightly faster plasma collisions than 
discussed here (0.5c and 0.9c instead of ~ 0.4c). 
In both cases the strong (n e /u! pe ~ 1) large-scale 
magnetic field is assumed to be obliquely oriented 
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relative to the flow axis, which provides for more ef- 
fective particle reflection in the perpendicular com- 
ponent of the magnetic field at the shock tran- 
sitions. Being interested in ion dynamics, Gar- 
gate & Spitkovsky (2012) present hybrid simula- 
tions and a PIC simulation of a parallel shock with 
Ma = 4.7, about 4.5 times the magnetization as- 
sumed here. They find emergent ion acceleration 
very late in the PIC simulation (at t ~ 87fl~ ) 
and also during the subsequent evolution captured 
with their hybrid simulations. They also find that 
ion acceleration proceeds to higher energies if the 
Alfvenic Mach number is larger, although the in- 
jection efficiency decreases. In ion plasma times, 
our magnetized PIC simulation has similar dura- 
tion as that of Gargate & Spitkovsky (2012), but 
is nearly 4 times shorter in ion gyro-times on ac- 
count of the lower plasma magnetization assumed 
in our simulations. However, the low Alfvenic 
Mach number used in the PIC simulation in Gar- 
gate & Spitkovsky (2012) suppresses Weibel-type 
instabilities. They are also not modeled with their 
hybrid simulations for faster shocks (their A/a = 31 
run should be comparable with ours). Because we 
demonstate that electron dynamics are important 
for the evolution of the system, a direct compar- 
ison with these results is not possible. Neverthe- 
less, all spectral structures we find coincide with 
anisotropic distribution functions, suggesting an 
imperfect ion-beam relaxation rather than accel- 
eration per se. 

• We have also studied the effect of the assumed 
ion-to-electron mass ratio on our results. We find 
that although the main characteristics of the long- 
time evolution of the systems do not critically de- 
pend on the mass ratio within the range studied, 
the detailed mechanisms of reaching a steady state 
might be modified or even different if the assumed 
mass ratio is too low. The timescales for particle- 
energy equilibration and the efficiency of particle 
pre-accelaration processes might be overestimated 
in simulations that use a reduced ion-to-electron 
mass ratio. 

Computer resource limitations prohibit us to study the 
properties of the shock formation for close-to-realistic 
ion-to-electron mass ratios. Studies of the saturation of 
the linear phase of purely-transverse ion-driven filamen- 
tation instability indicate a very low efficiency of the con- 
version of the ion beam energy into magnetic-field energy 
for a realistic mass ratio (e.g., Ren et al. 2007; Wiersma 
& Achterberg 2004), suggesting that other instabilities 
must mediate the shock formation in astrophysical set- 
tings. However, the applicability of these estimates to 



shocks might be limited, e.g., because they do not con- 
sider a nonlinear evolution of the filamentation instabil- 
ity and also not a steadily driven instability. In fact, re- 
cent numerical simulations of relativistic collisions of the 
electron-ion plasmas demonstrate the formation of the 
filamentation-instability mediated shocks, whose proper- 
ties do not significantly change with ion-to-electron mass 
ratio between 16 and 1000 (Spitkovsky 2008). Since the 
characteristics of the magnetic-field generation in our 
simulations are analogous to these observed in relativis- 
tic shocks, we similarly expect the properties of Weibel- 
instability-mediated nonrelativistic shocks will not be 
substantially modified for a realistic mass ratio. 

Concerning the energy equipartition between electron 
and ions, as well as the absence of significant non-thermal 
tails in the particle spectra, the question arises to what 
degree residual 2-body collisions in the simulation pro- 
vide thermalization that in the highly-collisionless space 
plasma would not be relevant. In space environments 
the mean free path for collective plasma interactions is 
many orders of magnitude smaller than that of 2-body 
collision. In PIC simulations one uses few computational 
particles to represent very many real electrons or ions 
and thus introduces artificial collisionality. Would that 
impact, and possibly prevent, particle pre-acceleration in 
our simulations? 

Coulomb collisions occur with a minimum impact pa- 
rameter of about 1 grid-cell length on account of the 
extended charge distribution used to represent the parti- 
cle. There are always many ions and electrons within one 
grid cell, and therefore substantial shielding will occur. 
Tests indicate that 2-body Coulomb interactions are in- 
deed negligible, but random charge imbalances that arise 
from a low number of particle per cell are possibly not. In 
fact, they can provide artificial scattering on time scales 
similar to or shorter than the simulation time for electron 
velocities < 0.03 c. The post-shock electrons, in partic- 
ular those in the spectral tail, are always considerably 
faster than that, and therefore we conclude that in our 
simulations artificial scattering can not suppress the for- 
mation of non-thermal tails in the electron spectrum. 

The authors thank Mark Dieckmann for comments 
and useful discussions. The work of J.N. is supported 
by Narodowe Centrum Nauki as research project DEC- 
2011/01/B/ST9/03183. M.P. and V.W. acknowledge 
support through grant PO 1508/1-1 of the Deutsche 
Forschungsgemeinschaft. The work of A.B. is supported 
by the project ENE2009-09276 of the Spanish Ministe- 
rio de Educacion y Ciencia. Simulations have been per- 
formed at the Pleiades facility at the NASA Advanced 
Supercomputing (NAS). 



APPENDIX 

SIMULATIONS WITH LOWER ION-ELECTRON MASS RATIO 

The present study is based on simulations that use a reduced ion-electron mass ratio, uii/m e = 50. To ascertain the 
effect of this parameter choice on our results, two additional large-scale numerical experiments with mi/m e = 20 have 
been performed (run U2 and M2; see Table 1). 

We find that the main characteristics of the long-time evolution of the systems do not critically depend on the ion- 
electron mass ratio within the range studied. In both the unmagnetized and magnetized case, a system of forward and 
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reverse shocks, separated by a CD, is formed. The late-stage structure and nature of the electromagnetic turbulence 
is analogous to that observed for runs with mi/m e — 50. Also, the processes of particle-energy dissipation and re- 
distribution proceed qualitatively in the same fashion. However, the detailed mechanisms of reaching a steady state 
might be modified or even different if the assumed mass ratio is too low. In the following, we only detail the most 
important modifications of the system characteristics introduced by a lower mass ratio. 



t=38Q3cd pI ~ 1 t = 3803oj pi ~ 1 




2500 3000 3500 2500 3000 3500 

x/A si x/X si 



Fig. 16. — Structure of the collision region at time t = 3803 uj^ 1 for (left) unmagnetized run U2 and (right) magnetized run M2 (see 
Fig. 3). In Figs. 16a and 16g, results for run Ul/Ml are plotted for t = 3803 ui~? with green lines (dense plasma) and blue lines (tenuous 
plasma) for comparison. The horizontal dashed lines mark the hydrodynamic compression level for the forward and reverse shocks (see 
Figs. 3 and 8). Figs. 16b and 16h overplot the magnetic and electric energy-density profiles for runs Ul and Ml, respectively. The mean 
magnetic field with energy density of 7.9 X lO - ^ is subtracted in the magnetic energy-density profile in Fig. 16h. 

Figure 16 presents the structure of the collision region at the end of run U2 (left) and M2 (right) at time t = 3803 w" 1 
(to be compared with Figs. 3 and 9). Overlaid on the average particle densities and electromagnetic-field profiles in 
Figures 16(a-b) and 16(g-h) are the profiles obtained for runs Ul and Ml. The horizontal axis is normalized to the 
skin depth of dense ions, which is X S i ~ 36.2 for mi/m e = 20. The particle spectra, measured in the CD rest frame, 
in the downstream region of the forward shock are shown in Figure 17 for run U2. Figure 18 displays spectra in the 
downstream region of the reverse (a) the forward (b) shock for run M2. 
Unmagnetized simulation: 

• The initial collision of the tenuous- and dense-plasma slabs leads to magnetic-field generation via Weibcl-likc 
instabilities, whose level is similar in run Ul and U2. However, due to their lower mass, the ions are more easily 
deflected in the magnetic field. Consequently, the ion-isotropization rate and the rate of kinetic-energy exchange 
with electrons are faster. Thus, in run U2 the CD acquires its steady-state velocity earlier, which accounts for a 
shift of the average profiles for runs Ul and U2 in Figures 16a and 16b. A substantial part of incoming beam 
particles is also reflected towards their respective upstream. In result, the contribution of dense ions to the 
forward-shock compression is significantly smaller in run U2 (compare solid red and green lines to solid black 
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tenuous electrons 




Fig. 17. — Normalized kinetic-energy spectra of tenuous electrons (black solid line), tenuous ions (red solid line), and the dense ions (red 
dashed line) downstream of the forward shock at x/\ s i m 3100 and time t = 3803 u>~^ for the unmagnetized run U2, in the CD rest frame. 
The ion-electron mass ratio is mi/m e = 20, and in simulation units E^ n = 0.25 (7 — 1) for electrons and -E^in = 5 (7 — 1) for ions. 




0.01 0.10 1.C 0.1 1.0 



L kin L kin 
Fig. 18. — (a) Kinetic-energy spectra of dense electrons (black line) and dense ions (red line) in the reverse-shock downstream region 
at x/\ s i ~ 2850 — 2890 and time t = 3803 uj^ 1 for the magnetized run M2. The dotted line shows a Maxwcllian fit to the electron 
spectrum, (b) Kinetic-energy spectra of tenuous electrons (black line) and tenuous ions (red line) in the forward-shock downstream region 
at x/X s i R2 3000 and time t = 3803^"^ for the magnetized run M2. The dotted line shows a relativistic Maxwellian for comparison. The 
spectra are shown in CD rest frame and are normalized and expressed in simulation units, in which m e c 2 = 0.25. 



and blue lines in Fig. 16a). However, also in this case the tenuous-ion filaments do not fully disintegrate in the 
shock downstream. Note a slight difference in the electron-compression profile right from the CD (dotted black 
and blue line in Fig. 16a) for the two mass ratios. 

• The forward-shock electromagnetic-turbulence structure in run U2 is analogous to that of run Ul. The role 
of the dense-ion beam in magnetic-field generation via Weibel-like instabilities and in the amplification of the 
parallel charge-density perturbations is in run U2 shared more equally between the beam of dense ions and the 
returning tenuous ions (compare Fig. 16e-f and 3e-f). The magnetic-turbulence level in the filamentation region 
of the shock precursor (3200 < x/\ s i < 3500) is considerably lower for m,/m e = 20, which results from the 
higher beam temperatures that cause a serious reduction of the growth rate of filamentation-like modes. 

• The tenuous-ion beam that permeates the dense-plasma region left from the CD is much warmer in run U2 
(compare Fig. 16f and 3f). This causes a decrease in the growth rate of the oblique filamentation modes and 
results in much weaker turbulent precursor to the reverse shock in run U2 (Fig. 16b). However, an increased 
efficiency of ion reflection in the CD region leads to an effective decoupling of the dense-ion beam from the 
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tenuous plasma; the numerous returning dense ions efficiently amplify the magnetic turbulence in the reverse- 
shock transition. The detailed structure of the reverse shock in run U2 is thus considerably modified with respect 
to that observed in run Ul (Fig. 16a-b for x/\ S i < 3010). However, also for mass ratio rrii/m e = 20, the energy 
density of the magnetic field left of the CD is more than an order of magnitude smaller than in the forward-shock 
region and the reverse-shock compression does not reach the hydrodynamic limit. 

• Particle spectra downstream of the forward shock (calculated in a slice of width ~ 20A S i, centered at x/X S i « 
3100, i.e., similar location as in run Ul) are shown in Figure 17. Whereas the electron and dense-ion spectra 
display characteristics analogous to those described for run Ul, the tenuous-ion distribution resembles that 
observed downstream of the forward shock in the magnetized simulation, Ml. This includes a plateau section in 
the momentum distribution, which only slightly deviates from isotropy in run U2. The contribution of tenuous 
ions to the total kinetic energy is, at a level of 38%, the same as in run Ul. The dense ions instead carry 17% 
only; the difference of 8% compared to run Ul is shifted to the electrons which now account for 44% of the 
total kinetic energy. These characteristics suggest that the timescales for particle-energy equilibration might be 
overestimated in lower-mass-ratio simulations. 

A slight dependence of the average density profiles on the ion-electron mass ratio has been reported in Kato & Takabc 
(2008). In our simulations, we see a similar effect in the forward-shock profile. Even more significant modifications 
are observed in the region to the left from the CD. However, the reverse-shock formation has not completed yet by 
the end of runs Ul and U2. 
Magnetized simulation: 

• A reduced mass ratio in run M2 results in larger Alfven velocities, which for the tenuous and dense plasma now 
read, «a,r = 2.9- 10 _2 c and «a,l = 8.8- 10 _3 c, respectively. Thus the simulation probes a slightly lower Alfvenic- 
Mach-number limit (but still Ma ^ 1), with Mach numbers Mj^i — 13.5 and M^ T ~ 14.4 for the forward and 
reverse shock, respectively. 

• The characteristic differences introduced by a lower ion mass, and hence larger rate of ion isotropization, that 
we described above for unmagnetized simulations, are also observed in the magnetized runs Ml and M2. This 
suggests that it is the ion-beam temperature, and not the slight change in the magnetic-field strength the ions 
experience through an increase of their cyclotron frequency, that is responsible for the observed differences in 
the development of Weibel-like instabilities in both runs. The level of magnetic turbulence in the forward-shock 
precursor in run Ml considerably exceeds that in run M2 (Fig. 16h). This is because in the magnetized case the 
turbulence is generated exclusively via filamentation-like instabilities, whose growth rates falls off with plasma 
temperature. Note, that the magnetic field decays faster downstream of the forward shock in run M2. 

• The effects of the warm beams on the system evolution are even more pronounced in the region left from 
the CD. Instead of the purely-transverse filamentation-like mode excited in the dense-plasma region by an 
energetic tenuous-ion beam in run Ml (see Sec. 3.2.1), a much weaker oblique mode is generated in simulation 
M2. Consequently, the current filaments are weak and their merging does not lead to the strongly-nonlinear 
"magnetic bubble" effects that in run Ml were largely responsible for the efficient decoupling of the colliding ion 
beams. Nevertheless, also in this case the system can form the reverse-shock transition at an early stage. The 
shock formation is now mediated by a strong filamentation-like instability due to the large number of dense ions 
reflected in the collision region. 

• As in run Ml, the formation of the reverse-shock transition at late stage in run M2 involves the amplification and 
compression of a parallel magnetic Bell-like mode. These processes are responsible for the turbulent-field structure 
at x/X S i < 2900 in Figure 16h. The parallel left-hand circularly polarized mode has properties similar to these 
observed in run Ml: the wavelength A ~ 20A S i, frequency u ~ 0.43 f2j, and the phase velocity v p h — 1.4i>a,l- 
Note, however, that the reverse-shock compression in run M2 does not reach the hydrodynamic limit by the end 
of the simulation. 

• Particle spectra shown in Figure 18 display characteristics similar to the distributions derived for run Ml. 
However, as is the case of run U2, the energy-equilibration processes proceed with a slightly faster rate for a 
lower mass ratio. The kinetic-energy content of the dense electrons downstream of the reverse shock is in run 
M2 the same within 0.3% as that of the dense ions, and the electrons are isotropic and thermal (compare a 
Maxwellian fit to the electron spectrum in Fig. 18a). Downstream of the forward shock, the tenuous electrons 
carry 45% of the total kinetic energy. An analysis of the time evolution of the electron spectra shows an indication 
of supra-thermal tails, that were not observed in run Ml. Therefore, although the electron spectra in run M2 
may be influenced by low statistics at high energies, this result demonstrates that the efficiency of particle 
injection/pre-accelaration processes might be overestimated in simulations that use a reduced ion-to-electron 
mass ratio. 
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TABLE 2 

Parameters a n and f) n of the dielectric 

TENSOR FOR THE SETUPS OF FIGURES 6 AND 12; 
RUN Ml. 



Figure 03 04 a 5 f}\ 



6 0.1 0.01 0.09 0.0354 0.354 

12 0.05 0.005 0.045 0.0354 0.3 



DIELECTRIC TENSOR ELEMENTS 
The elements of the dielectric tensor (1) read, 

1 + R (a 4 + a 5 + a 3 R)fi 2 
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In these expressions, Z x — k x VL/u> pet L is the reduced flow-aligned wave vector, R is the electron-to-ion mass ratio, and 
£Ib — Sle/wp e ^L the reduced electron cyclotron frequency. The parallel magnetic mode appears in the reverse-shock 
precursor in the late stage of our magnetized simulations Ml and M2 (Fig. 12). However, the mode can also be 
observed in the linear spectrum of the early-stage turbulence (Fig. 6). We provide parameters a n and j3 n for both 
setups in Table 2. 



REFERENCES 



Hoshino, M. 2009, ApJ, 690, 244 

Hoshino, M. 2010, Physical Review Letters, 104, 



Amano, T.. 
Amano, T.. 
181102 

Bell, A.R. 2007, MNRAS, 353, 550 

Blandford, R., & Eichler, D. 1987, Phys. Rep., 154, 1 

Bret, A. 2007, Computer Physics Communications, 176, 362 

Bret, A. 2009, ApJ, 699, 990 

Bret^ A., Firp'o, M.-C, & Deutsch, C. 2005, Phys. Rev. Lett., 94, 
115002 

Bret, A., Dieckmann, M. E., & Deutsch, C. 2006, Phys. Plasmas, 
13, 082109 

Bret, A., Gremillet, L., & Benisti, D. 2010, Phys. Rev. E, 81, 
036402 

Buneman, O. 1993, in Computer Space Plasma Physics: 
Simulation Techniques and Software, Eds.: Matsumoto & 
Omura, Tokyo: Terra, p. 67 
Burgess, D., & Scholer, M. 2007, Physics of Plasmas, 14, 012108 
Cargill, P. J., & Papadopoulos. K. 1988, ApJ, 329, L29 
Chang, P., Spitkovsky, A., & Arons, J. 2008, ApJ, 674, 378 
Classen, H.-T., & Mann, G. 1997, A&A, 322, 696 
Dieckmann, M. E., Kourakis, I., Borghesi, M., & Rowlands, G. 

2009, Phys. Plasmas, 16, 074502 

Dieckmann, M. E., Murphy, G. C, Meli, A., & Drury, L. O. C. 

2010, A&A, 509, A89 

Fiore, M., Silva, L. O., Ren, C, Tzoufras, M. A., & Mori, W. B. 

2006, MNRAS, 372, 1851 
Fonseca, R. A., Martins, S. F., Silva, L. O., Tonge, J. W., Tsung, 
F. S., & Mori, W. B. 2008, Plasma Phys. and Control Fusion, 
50, 124034 

Freidberg, J. P. 1971, Phys. Rev. Lett., 27, 



Forslund, D. W., 

1189 
Forslund, D. W., 
Frederiksen, J. T 



Shonk, C. R. 1970, Phys. Rev. Lett., 25, 1699 



Hededal, C. B., Haugboelle, T., & Nordlund, 
A. 2004, ApJ, 608, L13 
Gargate, L., & Spitkovsky, A. 2012, ApJ, 744, 67 
Gary, S. P. 1984, J. Geophys. Res., 89, 179 
Goldman, M. V. 1984, Rev. Mod. Phys., 56, 709 
Greenwood, A. D., Cartwright, K. L., Luginslund, J. W., & Baca, 
E.A. 2004, J. Comp. Phys., 201, 665 



Frederiksen, J. T., & Nordlund, 



Hededal, C. B., Haugboelle, T.. 

A. 2004, ApJ, 617, L107 
Kazimura, Y., Califano, F., Sakai, J. -I., et al. 1998b, Journal of 

the Physical Society of Japan. 67, 1079 
Kazimura, Y., Sakai, J. I., Neubert, T., & Bulanov, S. V. 1998a, 

ApJ, 498, L183 
Kato, T. N., & Takabe, H. 2008, ApJ, 681, L93 
Kato, T. N., & Takabe, H. 2010, ApJ, 721, 828 
Krauss-Varban, D., & Omidi, N. 1991, J. Geophys. Res., 96, 17715 
Lembege, B., Giacalone, J., Scholer, M., et al. 2004, 

Space Sci. Rev., 110, 161 
Lembege, B., Savoini, P., Hellinger, P., & Travmcek, P. M. 2009, 

J. Geophys. Res., 114, 3217 
Matsumoto, Y., Amano, T., & Hoshino, M. 2012, arXiv:1204.6312 
McClements, K. G., & Fletcher, L. 2009, ApJ, 693, 1494 
Medvedev, M. V. 2006, ApJ, 651, L9 
Medvedev, M. V., & Loeb, A. 1999, ApJ, 526, 697 
Medvedev, M. V., Fiore, M., Fonseca, R. A., Silva, L. O., & Mori, 

W. B. 2005, ApJ, 618, L75 
Murphy, G. C, Dieckmann, M. E., Bret, A., & Drury, L. O'C. 



Drury, L. 
& Stroman, T. 2010, ApJ, 709, 



z Nishikawa, K.-I. 2008, ApJ, 
Richardson, G., Preece, R., Sol, 



Fishman, G. J. 



2010, A&A, 524, 
Niemiec, J., Pohl, M., Bret, A. 
1148 

Niemiec, J., Pohl, M., Stroman, 

684, 1174 
Nishikawa, K.-I., Hardee, P. E. 

H., & Fishman, G. J. 2003, ApJ, 595, 555 
Nishikawa, K.-I., Hardee, P. E., Hededal, C. B.. 

2003, ApJ, 642, 1267 
Nishikawa, K.-I., Niemiec, J., Hardee, P. E.. 

Sol, H., Mizuno, Y., Zhang, B., Pohl, M„ 

Hartmann, D. H. 2009, ApJ, 698, L10 
Papadopoulos, K. 1988, Ap&SS, 144, 535 
Polomarov, O., Kaganovich, I., & Shvets, G. 

Phys. Rev. Lett., 101, 175001 
Quest, K. B. 1988, J. Geophys. Res., 93, 9649 
Ren, C, Blackman, E. G., & Fong, W. 2007, Phys. Plasmas, 17, 

012901 



Medvedev, M. 
Oka, M., & 



2008, 



V.. 



22 



Nicmiec et al. 



Reville, B., Kirk, J. G., Duffy, P. 2006, Plasma Physics and 

Controlled Fusion, 48, 1741 
Riquelme, M. A., & Spitkovsky, A. 2009, ApJ, 694, 626 
Riquelme, M. A., & Spitkovsky, A. 2011, ApJ, 733, 63 
Shvets, G., Polomarov, O., Khudik, V., Siemon, K., & 

Kaganovich, I. 2009, Phys. Plasmas, 16, 056303 
Silva, L. O., Fonseca, R. A., Tonge, J. W., Dawson, J. M., Mori, 

W. B., Medvedev, M. V. 2003, ApJ, 596, L121 
Spitkovsky, A. 2008, ApJ, 673, L39 
Stasicwicz, K., & Ekeberg, J. 2008, ApJ, 680, L153 
Stroman, T., Pohl, M., Niemiec, J., Bret, A. 2012, ApJ, 746, 24 
Stroman, T., Pohl, M., & Niemiec, J. 2009, ApJ, 706, 38 



Sugiyama, T. 2011, Phys. Plasmas, 18, 022302 

Tzoufras, M., Ren, C, Tsung, F. S., Tonge, J. W., Mori, W. B., 

Fiore, M., Fonseca, R. A., & Silva, L. O. 2006, 

Phys. Rev. Lett., 96, 105002 
Umeda, T., Yamao, M., & Yamazaki, R. 2009, ApJ, 695, 574 
Weibel, E. S. 1959, Phys. Rev. Lett., 2, 83 
Wiersma, J., & Achterbcrg, A. 2004, A&A, 428, 365 
Winske, D., Leroy, M.M. 1984, J. Geophys. Res., 89, 2673 
Winske, D., Thomas, V. A., Omidi, N., & Quest, K. B. 1990, 

J. Geophys. Res., 95, 18821 



